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ABSTRACT 



We evaluate the construction methodology of an all-sky catalogue of galaxy clusters detected through the Sunyaev-Zel'dovich (SZ) 
^ ■ effect. We perform an extensive comparison of twelve algorithms applied to the same detailed simulations of the millimeter and 

\^ \ submillimeter sky based on a Planck-like case. We present the results of this "SZ Challenge" in terms of catalogue completeness, 

purity, astrometric and photometric reconstruction. Our results provide a comparison of a representative sample of SZ detection 
■ algorithms and highlight important issues in their application. In our study case, we show that the exact expected number of clusters 

remains uncertain (about a thousand cluster candidates at \b\ > 20 deg with 90% purity) and that it depends on the SZ model and on 
_ _ the detailed sky simulations, and on algorithmic implementation of the detection methods. We also estimate the astrometric precision 

' , of the cluster candidates which is found of the order of ~ 2 arcmins on average, and the photometric uncertainty of order ~ 30%, 

depending on flux. 
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' 1. Introduction larger volumes to higher redshifts with vastly superior statistics 

^ , and control of systematics. 

The Planck SZ cluster catalogue will be one of the im- 

Galaxy cluster catalogues have played a long-standing, vital portant players in this context. Surveying the entire sky in 9 

role in cosmology, providing important information on top- millimeter/sub millimeter bands with ~ 5 - 10 arcmin resolu- 

ics r anging from cosmological par ameters to galaxy forma- tion over the channels most sensitive to the cosmic microwave 

tion (iRosati et alJ l2002t IVoitI l2005h . In particular, recent X- background (CMB) anisotropies, the Planck satellite will find 

ray cluster catalogs have proved valuable in establishing the large numbers of clusters through the Sunyaev-Zel'dovich (SZ) 

standard cosmo logic al model (e.g., Schuecker et al. (2003), effect (Sunyaev & Zeldovich 1970, 1972; Birkinshaw 1999|; 

IVikhlinin et aH (12(1091) 1. The science potential of large cluster ICarlstrom et alJ200 2). The advantages of this, much anticipated, 

surveys is strong: They are, for instance, considered one of the technique include efficient detection of distant clusters and se- 

central observational tools for illuminating the nature of dark lection based on an observable expected to correlate t ightly with 

energv (e.g.. the Dark Energv Task Force Report (lAlbrecht et alJ cluster mass (Bartlett 2002; da Silva et al. 2004a'; Motl et ajj 

[2006, 2009)). A suite of lai-ge cluster surveys planned over the 2005a; Nagai 2006; Bonaldi et al. 2007; Shaw et al. 2008). An 

coming years in the optical/IR, X-ray and millimeter bands will official mission deliverable, the Planck SZ catalogue will be the 

greatly extend the reach of cluster science by probing much first all-sky cluster catalogue since the workhorse catalogs from 
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the ROSAT All-Sky Survey (RASS) (lTruemnedll992h : in other 
words, Planck will be the first all-sky cluster survey since the 
early 1990s! 

Within the Planck Consortium, a considerable effort has been 
conducted for the scientific evaluation of the cluster catalogue 
construction methodology. As part of this evaluation effort, we 
completed an extensive comparison of twelve algorithms applied 
to detailed simulations of Planck data based on the Planck Sky 
Model. This study was dubbed "The SZ Challenge" and was car- 
ried out in two steps using different SZ cluster models and cos- 
mologies; these are referred to as Versions 1 and 2 and more fully 
explained below. We report the findings of these initial studies in 
terms of catalogue completeness and purity, as well as astromet- 
ric and photometric accuracy and precision. 

The article is organized as follows: In the next section we 
detail our sky simulations, including a brief description of the ir 
basis, the Planck Sky Model (PSM, Del abrouille et all (120121) '). 
The following section then introduces the different catalogue 
construction methodologies employed, before moving on to a 
presentation of each of the twelve algorithms in the study. We 
present the results of the challenge in Section |4] followed by 
a comparative study of algorithmic performance. We conclude 
with a discussion of both the limitations of this study and future 
directions. 



2. Simulations 

2.1. The Planck Sky Model 

One of the important differences between the present work and 
previous studies of Planck SZ capabilities is the detailed and 
rather sophisticated simulation of milUmeter and submillime- 
ter sky emission used here. Our sky simulations are based on 
an early development ve rsion of the Planck Sky Model (PSM, 
iDelabrouille et al.l (120121) ). a flexible software package devel- 
oped within the Planck collaboration for making predictions, 
simulations and constrained realisations of the microwave sky. 
The simulations used for this challenge are not polarised (only 
temperature maps are useful for detecting clusters using the ther- 
mal SZ effect). 

The CMB sky is based on the best fit angular power spectrum 
model of WMAP The CMB realisation is not constrained by 
actual observed CMB multipoles, in contrast to the simulations 
used by Leach et al. (2008). 

Diffuse Galactic emission is described by a four compo- 
nent model of the interstellar medium comprising free-free, 
synchrotron, thermal dust and spin ning dust and is based on 
iMiville-Deschenes et al] ( |2008|) (see iMiville-D eschenes' ( |2009|) 
for a review). The predictions rely on a number of sky tem- 
plates with different angular resolutions. In order to simulate 
the sky at Planck resolution we have added small-scale fluctu- 
ations to some of the templates to increase the fluctuation level 
as a function of the local brightness and therefore reproduce 
the non-Gaussian and non-uniform properties of the interstellar 
emission. The procedure used t o add small scales is presented in 
iMiville-Deschgnes et al] (l2007h . 

Fre e-free emission is based on the model of lDickinson et alj 
(I2003h assuming an electronic temperature of 7000 K. The spa- 
tial structure of the emission is modeled using an Ho- tem- 
plate corrected for dust extinction. The Ha map is a combi- 
nation of the Southern H-Alpha Sky Survey Atlas (SHASSA, 
[Oaustad e t al. (2001)) and th e Wisconsin H- Alpha Mapper 
(WHAM. lHaffneret alJ dlool ). smoothed to obtain a uniform 
angular resolution of 1°. Dust extinction is inferred using the 



E(B - V) all-sky map of lSchlegel et"an (Il998h . As mentioned 
earlier, small scales were added in both templates to match the 
Planck resolution. The free-free emission law is constant over 
the sky, as it depends only on the electronic t emper ature, taken as 
a constant here (see however [Wakker et al.l (l2008l) for a descrip- 
tion of high-velocity clouds not detected by the WHAM survey 
and hence not included in our simulations). 

Synchrotron emission is based on an extrapolation of the 
408 MHz map of iHaslam etall ( Il982h from which an estimate 
of the free-free emission was removed. In any direction on the 
sky, the spectral emission law of the synchrotron is assumed to 
follow a power law, j'*'"'^ oc yP. We use a pixel-dependent spec- 
tral index (3 derived from the ratio of the 408 MHz map and the 
estimate of the s ynchrotron emission at 23 GHz in the WMAP 
data obtained bvMivil le-Deschenes et alj ( l2008l) . 

The thermal emissi on from i nterst ellar dust is estimated us- 
ing model 7 of Finkbei ner et a (119991) . This model, fitted to the 
FIRAS data (7° resolution), assumes that each line of sight can 
be modeled by the sum of the emission from two dust popula- 
tions, one cold and one hot. Each grain population is in thermal 
equilibrium with the radiation field and thus has a grey-body 
spectrum, so that the total dust emission is modelled as 

2 

I^^Yjf^^'BviTd (1) 

where B^iTi) is the Planck function at temperature T, . In model 
7 the emissivity indices are fi\ - 1.5, - 2.6, and f\ - 0.0309 
and /2 - 0.9691. Once these values are fixed, the dust temper- 
ature of the two components is determined using only the ratio 
of the observations at 100 jum and 240 ^m. For this purpose, 
we use the 100/240 /im map ratio published bv lFinkbeiner et aP 
('1999'). Knowing the temperature and y6 of each dust component 
at a given position on the sky, we use the 100 jum brightness at 
that position to scale the emission at any frequency using Eq. ([T]i. 

Spinning dust emission uses as a template the dust extinction 
map E{B - V), and us es an emission law uniform on the sky, 
based on the model of iDraine & LazarianI (Il998h . assuming a 
Warm Neutral Medium (WNM). 

We emphasise that the emission laws of both synchrotron 
and dust vary across the sky. The spectral index of free-free and 
the emission law of spinning dust, however, are taken as uniform 
on the sky. 

Point sources are modeled with two main categories: radio 
and infrared. In the present simulation, none of two is corre- 
lated with the SZ signal. Simulated radio sources are based on 
the NVSS or SUMSS and GB6 or PMN catalogues. Measured 
fluxes at 1 and/or 4.85 GHz are extrapolated to 20 GHz using 
their measured SED when observed at two frequencies. Sources 
for which a flux measurement is available at a single frequency 
have been randomly assigned to either the steep- or to the flat- 
spectrum class i n the p roportions observationally determined by 
iFomalont et aP (119911) for various flux intervals, and assigned 
a spectral index randomly drawn from the corresponding dis- 
tribution. Source counts at 5 and 20 GHz obtained in this way 
are compared, for consistency, with observed counts, with the 
model by Toffolatti et al. ( 1998), and with an updated version of 
the model by de Zotti et al.l (120051) . allowing for a high-redshift 
decline of the space density of both flat-spectrum quasars (FSQs) 
and steep-spectrum sources (not only for FSQs as in the origi- 
nal model). Further extrapolation at Planck frequencies has been 
made allowing a change in SED above 20 GHz, assuming again 
a distribution in flat and steep populations. For each of these 
two populations, the spectral index is randomly drawn within a 
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set of values compatible with the typical average and dispersion. 
Infrared sources are based on the IRAS catalogue, and modeled 
as dusty galaxies (Serjeant & Harrison 2005). IRAS coverage 
gaps were filled by randomly adding sources with a flux distri- 
bution consistent with the mean counts. Fainter sources were as- 
sumed to be mostly sub-millimeter bright galaxies, such as those 
detected by SCUBA surveys. These were modelled following 
[Ora nato et al. (2004) and assumed to be strongly clustered, with 
a comoving clustering radius ro - 8 h"' Mpc. Since such sources 
have a very high areal density, they are not simulated individu- 
ally but make up the sub-mm background. 

The SZ component is described in detail in the following 
section. 

Component maps are produced at all Planck central frequen- 
cies. They are then co-added and smoothed with Gaussian beams 
as indicated in Table[T] extracted from the Planck Blue book. We 
thus obtain a total of nine monochromatic sky maps. 

Finally, inhomogeneous noise is simulated according to the 
pixel hit count corresponding to a nominal 14-month mi ssiorfl 
using the Level-S simulation tool (iReinecke et al.l l2b06h . The 
RMS noise level in the simulated maps is given in Table[T]from 
the Planck Blue book. It is worth noting that the in-flight perfor- 
mances of Planck-HFI as reported by the Planck collaboration 
dPlanck Collaboration et al .11201 lal) are better than the require- 
ments. 



2.2. The SZ component 

We simulate the SZ component using a semi-analytic approach 
based on an analytic mass function dN{M,z)/dMdz. After se- 
lecting cosmological parameters (h, Q.^, Q.^, erg), the cluster 
distribution in the mass-redshift plane (M, z) is drawn from a 
Poisson law whose mean is given by the mass function. Clusters 
are spanning the mass range 5 x 10'-' Mq < Mvh- < 5 x IO'^Mq 
and the redshift range < z < 4. Cluster Galactic coordinates 
(I, b) are then uniformly drawn on the sphere. We compute the 
SZ signal for each dark matter halo following two different mod- 
els, producing two simulations (vl and v2) with different sets of 
cosmological parameters and mass functions and SZ signal^ 



2.2.1. SZ Challenge vl 

For the first version of the SZ Challenge, we used h - 0.7, 
- .3, - 0.7, cTg - 0.8 5 and the Sheth-Tormen mass 
function ISheth & Tormei] (Il999h . We assume that the clusters 
are isothermal and that the electron density profile is given by 
the y6-model, with fi - 2/3, and core radius scaling as M'^-'. We 
truncate the model at the virial radius, ryir, and choose the core 
radius - rvir/lO. The virial radius here is defined according to 
the spherical collapse model. The temperat ure of each cluster is 
derived using a mass-temperature given in IColafrancesco et aP 
withTjs - 7.75keV, consistent with the simulations of 
lEkeetal.l(ll998h. For more details on this model, we refer the 
reader to § 5 of Ide Zotti et al.1 (l2005h . 



2.2.2. SZ Challenge v2 

The second version of the S Z Challenge was pr oduced using 
a WMAP5 only cosmology (iDunklev et all 120091) {h = 0.719, 
= 0.256, i^4_^=„0J44^8 = 0.798). We used the Jenkins 
mass function (fjenkins et al.ll200l1) . The SZ emission is mod- 
eled using the universal pre ssure profile d erived from the X-ray 
REXCESS cluster sample (lAmaud et al.l 12010) which predicts 
profile and normalisation of SZ clusters given their mass and 
redshift. The profile is well fitted by a generalized NEW profile 
that is much steeper than the y6-profile in the outskirts. Moreover, 
for a given mass, the normalisation of the SZ flux is ~ 15% 
lower than the normalisation of SZ Challenge vl. This profile 
was used as the baseline p rofile i n the SZ early results from 
Planck (Planck Collaboration et al.ll201 lb f.e.c. j) . 

Neither of the two sets of simulations (vl and v2) contains 
point sources within clusters. The eff'ect of contamination by ra- 
dio or infrared point sources in clusters was therefore not stud- 
ied her^l We neither include relativistic electronic populations 
within clusters. As for point sources, this effect was not studied 
herdl 



3. Methods and Algorithms 

The SZ Challenge was run as a blind test by providing the simu- 
lated sky maps. Participating teams, ten, were then asked to run 
their algorithms, twelve in total on the simulated data and supply 
a cluster catalogue with 

1 . (a, 6): cluster sky coordinates 

2. Yrec'- recovered total SZ flux, in terms of the integrated 
Compton-F parameter 

3. AYrec'- estimated flux error, i.e., the method's internal esti- 
mate of flux error 

4. Oiec'- recovered cluster angular size, in terms of equivalent 
virial radius 

5. AOrec'- estimated size error (internal error). 

The different methods were divided into two classes: direct 
methods that produce a cluster catalogue applying filters directly 
to a set of frequency maps, and indirect methods that first extract 
a thermal SZ map and then apply source finding algorithms. 

In this classification, the 12 algorithms studied were: 

- Direct methods: 

• MMFl: Matched Multi-Filters (MMF) as implemented 
by D. Harrison. 

• MMF2: MMF as described inlHe rranz etall (|2002|) . 

• MMF3: MMF as described in Meli n et al.1 (12006 ' 

• MMF4: MMF as described in Scha fer et all (12006 

• P wS: Bayesian meth od PowellSnakes as described 
in lCarvalho et al.l(l2009h . 

- Indirect methods: 



' Note that since launch in May 2009, the observed Helium consump- 
tion for the Planck-HFI dilution cooler indicates that the instrument can 
operate for about 30 months. A mission extension has been approved 
by ESA accordingly. 

^ We assume spherical symmetry for the individual SZ clusters and 
do not take into account any scatter in the distribution of pressure pro- 
files. 



^ The eff^ect of radio sources (y < 217GHz) is to reduce the ob- 
served SZ signal at a given frequency while the effect of infrared sources 
(y > 217GHz) is to increase it. However, the extraction algorithms be- 
ing multifrequency, their sensitivity to point sources is expected to be 
weaker than for single frequency extractions because of the different 
spectral dependence of point sources and SZ clusters. 

* The effect has been very recently studied within the Planck 
Collaboration: assuming a non relativistic spectrum for extracting clus- 
ters biases the flux low by about 10% in direction of massive (Myi,. > 
lO'^Mo) clusters (Planck Collaboration et al. 2013, in preparation). 
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Table 1. Characteristics of instrumental values taken from Planck's Blue Book for a 14 month nominal mission, o-pixei refers to the 
standard deviation of the 1.7' (Healpix nside=2048) pixel noise maps at the considered frequency. 



Channel 


30 GHz 


44 GHz 


70 GHz 


100 GHz 


143 GHz 


217 GHz 


353 GHz 


545 GHz 


857 GHz 


FWHM [arcmin] 


33 


24 


14 


10 


7.1 


5 


5 


5 


5 


0^ pixel [mA'cMs] 


0.131 


0.130 


0.126 


0.057 


0.029 


0.046 


0.137 


1.241 


56.639 



BNP; Ba yesian No n-parametric method as de- 
scribed in Diego et al. (2002), followed by SExtractor 



(iBertin & ArnoutSnl996.) to detect clusters and perform 
photometry. 

• ILCl: All-sky Internal Linear combination (ILC) on 
needlet coefficients^ ( similar to th e method used for 
CMB extraction in [Pelabrouille et al.1 (I2009D ) to get an 
SZ map, followed by Matched Filters on patches to ex- 
tract clusters and per form photometry. 

• ILC2: Same SZ map 'Pelabrouille et al.1 (I200 9'). but fol- 
lowed by SExtractor on patches instead of a matched fil- 
ter to extract the clusters. Photometry or flux measure- 
ment is however done using Matched Filters at the posi- 
tion of the detected clusters. 

• ILC3 developed by G. Chon and R. Kneissl: ILC in real 
space and filtering in harmonic space to obtain an SZ 
map, followed by fitting a cluster model. 

• ILC4 developed by J.-B. Melin: ILC on patches in 
Fourier space to obtain SZ maps, followed by SExtractor 
to detect clusters and Matched Multifilters to perform 
photometry. 

• ILC5 developed by D. Yvon: ILC on patches in wavelet 
space to obtain SZ maps, followed by SExtractor to de- 
tect clusters and perform photometry. 

• GMCA: Generalized Morphologi c al Component 
Analysis as described in iBobin et akl (l2008l) . followed 
by wavelet filtering and SExtractor to extract clusters 
and Matched Multifilters to perform photometry. 

All of the algorithms make use of the known frequency spec- 
trum of the SZ signal; attempts to detect clusters without this 
prior knowledge perform significantly worse. A summary of the 
characteristics of the codes as well as their treatment of the point 
sources, foreground removal and masking is given in Table |2] 
Further details about the algorithms are given in the following 
subsections. 



3.1. Matched Multifilter Methods 

The multifrequency matched filter (MMF) enhances the contrast 
(signal-to-noise) of objects of known shape and known spectral 
emission law over a set of observations containing correlated 
contamination signals. It offers a practical way of extracting a 
SZ clusters using multifrequency maps. The method makes use 
of the universal thermal SZ effect frequency dependance (as- 
suming electrons in clusters are non-relativistic), and adopts a 
spatial (angular profile) template. The filter rejects foregrounds 
using a linear combination of maps (which requires an estimate 
of the statistics of the contamination) and uses spatial filtering 
to suppress both foregrounds and noise (making use of prior 
knowledge of the cluster profile). In all cases discussed here, the 
adopted template is identical to the simulated cluster profiles, ex- 



^ needlet coefficients are the equivalent of Fourier coefficients in the 
adopted spherical wavelet domain 



cept for MMF 2 on SZ Challenge v2 . The MMF has been stu died 
extensively bv lHerranz et al.1 (l2002h and lMehn et~an (l2006h . 

Three of the MMF methods tested here work with projected 
flat patches of the sky, and one method works directly on the 
pixelised sphere. 

In the first case, the full-sky frequency maps are projected 
onto an atlas of overlapping square flat regions. The filtering is 
then implemented on sets of small patches comprising one patch 
for each frequency channel. For each such region, the nine fre- 
quency maps are processed with the MMF. A simple threshold- 
ing detection algorithm is used to find the clusters and produce 
local catalogues. The MMF is applied with varying cluster sizes 
to find the best detection for each cluster. This provides an esti- 
mate of the angular size in addition to the central Compton pa- 
rameter. Each algorithm explored its own, but similar, range of 
angular scales; MMF3, for example, runs from 6v - 2to 150 ar- 
cmins. The catalogues extracted from individual patches are then 
merged into a full-sky catalogue that contains the position of the 
clusters, their estimated central Compton parameter, the virial 
radius and an estimation of the error in the two later quantities. 
The integrated Compton parameter is derived from the value of 
the central Compton parameter and the radius of the cluster. 

In the following, we give relevant details specific to each im- 
plementation of the MMF. 

3.1.1. MMF1 

The performance of the MMFl algorithm is sensitive to the ac- 
curacy of the evaluation of the power spectra and cross-power 
spectra of the non-SZ component of the input maps. The de- 
tection is performed in two passes, the first detecting the high- 
est signal-to-noise SZ clusters, and the second detecting fainter 
clusters after the removal of the contribution of the brightest ones 
from the power spectra estimated on the maps. 

The merging of the catalogues from distinct patches is im- 
plemented with the option of discarding detections found in the 
smallest radius bin. These detections essentially correspond to 
spatial profiles indistinguishable from that of a point source. 
This option permits better control of the contamination by point 
sources, as a disproportionate fraction of the spurious detections 
occur in this bin ; despite their different spectral signature, point 
sources can occasionally pass through the filter. Using this op- 
tion of MMFl reduces the contamination at a given threshold 
(see Section 14.11 for definitions of this and other diagnostics of 
catalogue content) depending on the actual profile of the SZ clus- 
ters. 



3.1.2. MMF2 

The MMF2 a l gorith m follows closely the method described in 
iHerranz et al.l (l2002h . The method is simple and quite robust, 
although the performance depends on the model assumed for 
the radial profile of the clusters. For this work, a truncated mul- 
tiquadric profile similar to a y6-model has been used for SZ 
Challenge vl and v2. The profile is not a good match for the 
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Table 2. Summary of the algorithms compared in the SZ Challenge. The first column shows the name of the method. The second 
column indicates when the code is using a prior on the SZ cluster shape. The uperscript " indicate that the detection did not use a 
shape prior but that the computation of the SZ flux did. The third column gives the performance in terms of CPU hours needed to 
complete the analysis. The fourth and fifth column show whether the analysis was made using all-sky maps or projected patches, 
their area in square degrees and the number of patches. Methods ILCl and ILC2 work with full sky maps for producing an SZ map 
(by ILC on needlet coefficients) and then work with 504 small patches for cluster detection by matched filtering (ILCl) or using 
SExtractor (ILC2). The sixth and seventh columns provide information about any specific method used for subtracting point sources 
(PS) and Galactic foregrounds (FG). Note that both the MMF and the ILC methods have a built-in way for subtracting both point 
sources and diffuse foregrounds, by treating them as additional noise (of astrophysical origin) correlated across the channels. Note, 
also, that the study is made only on clusters at \b\ > 20 degrees Galactic latitude for all methods. The eighth column summarizes the 
main characteritics of each algorithm in terms of yield at 90% purity and photometric accuracy. 



Methiod 


Shape matching 


CPU 


Patches 


Number of 


PS subtraction 


FG subtraction 


Main characteristics 






Ximc (h) 


(size deg.) 


patches 


method 


method 




MMFl 


Yes 


50-60 


14.6° X 14.6° 


640 






Best yield among MMFs 

vJUUU. UIUJUJlilCLiy 


MMF2 


Yes 


31 


14.6° X 14.6° 


371 






Good yield 
Good photometry 


MMF3 


Yes 


5 


10° X 10° 


504 


mask lOo- PS 




Good yield 

vJUUU. UIUJUJlilCLiy 


MMF4 


Yes 




full sky 








Poor yield (see Sec. 13. 1.4b 
No photometry 


PwS 


Yes 


5.73 


14.6° X 14.6° 


2064 






Good yield 
Good photometry 


BNP 


No 


15 


10° X 10° 


512 


MHW 


Subtract 857 


Median yield 
Poor photometry 


ILCl 


Yes 


2-3 


see caption 


(504) 




Good yield 
Good photometry 


ILC2 


No 


2-3 


see caption 


(504) 




Best yield among ILCs 
Good photometry 


ILC3 


Yes 


24 


full sky 






Template fitting 


Poor yield 
Good photometry 


ILC4 


Yes" 


6 


10° X 10° 


504 


mask lOcr PS 




Good yield 
Good photometry 


ILC5 


No 


0.2 


11° X 11° 


461 


SExt. 




Poor yield (see Sec. 13.3.51 
Poor photometry 


GMCA 


No 


4 


10° X 10° 


504 






Median yield 
Good photometry 



simulated profile in SZ Challenge v2. This does not affect sig- 
nificantly the completeness and purity of the method (as shown 
in Section nil but the extracted flux is biased with respect to the 
input. The family of profiles used by the algorithm can however 
be adjusted. 



3.L3. MMF3 

The MMF3 SZ extraction algorithm is an all- sky extension of th e 
Matched Multifrequency Filter described in iMeUn et all (l2006h . 
It has been used for the production of the early SZ cluster sam- 
ple (Planck Collaboration et al. 2011b). In the version used for 
the SZ challenge, auto and cross power spectra used by the filter 
do not rely on prior assumptions about the noise, but are directly 
estimated from the data. They are thus adapted to the local in- 
strumental noise and astrophysical contamination. 



3.L4. MMF4 



The s pherical matched and scale adaptive filters dSchafer et al 
|2006|) are generalisations of the filters proposed bv lHerranz et al 



(1200 2') for spherical coordinates. Just like their counterparts they 
can be derived from an optimisation problem and maximise the 
signal to noise ratio while being linear in the signal (matched fil- 
ter) and being sensitive to the size of the object. The algorithm 



interfaces to the common HEALPIX package and treats the en- 
tire celestial sphere in one pass. 

The most important drawback is the strong Galactic contam- 
ination - the filter was not optimised to deal with a Galactic cut 
like many of the other algorithms, although it is in principle pos- 
sible to include that extension. The large noise contribution due 
to the Galaxy is the principal reason why the performance of the 
filter suffers in comparison to the approach of discarding a large 
fraction of the sky. 

3.2. Bayesian Methods 
3.2. L PowellSnakes 

PoweflSnakes (PwS) (iCarvalho et al.l (l2009l) ) is a fast multi- 
frequency Bayesian detection algorithm. It analyses flat sky 
patches using the ratio 



P = 



Pi(Hi\d) 



(2) 



where Hi is the detection hypothesis, 'There is a source' and Hq 
the null hypothesis 'Only background is present' CJaynes 2003i) . 
Applying Bayes theorem to the above formula one gets 



Pr(rf|j/i)Pr(//i) ^ Zi Pr(//i) 
Pr(rf|//o)Pr(//o) ZoPi-(//o)' 



(3) 
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where 

j £(Q)n(0)d^Q, (4) 

is the evidence, -£(©) is the hkehhood, 7r(0) is the prior and a 
vector representing the parameter set. 

An SZ parameterised template profile of the clusters 
s(X, A) = ^ T{a, X - X, y - Y), is assumed known and fairly 
representative of the majority of the clusters according to the 
resolution and S/N ratio of the instrumental setup, where t(. . .) 
is the general shape of the objects (beta or Arnaud et al. profile) 
and Uj a vector which contains the parameters controlling the 
geometry of one specific element (core/scale radius, parameters 
of the beta or Arnaud et al. profile). 

The algorithm may be operated on either 'Frequentist mode' 
where the detection step closely resembles a multi-frequency 
multi-scale 'Matched Filter' or 'Bayesian mode' where the pos- 
terior distributions are computed resorting to a simple 'Nested 
Sampling' algorithm CFeroz & Hobsoa,2008,) . 

The acceptance/rejection threshold may be defined either by 
using 'Decision theory' where the expected loss criterion is min- 
imised or by imposing a pre-defined contamination ratio. In the 
case of a loss criterion, the symmetric criterion — 'An unde- 
tected cluster is as bad as spurious cluster' — is used. 

3.2.2. BNP 

This method is described in detail in iDiego et alJ (l2002h . It is 
based on the maximization of the Bayesian probability of having 
an SZ cluster given the multifrequency data with no assumptions 
about the shape nor size of the clusters. The method devides the 
sky into multiple patches of about 100 sq. degrees and performs 
a basic cleaning of the Galactic components (by subtracting the 
properly weighted 857 GHz map from the channels of inter- 
est for SZ) and point sources (using a Mexican Hat Wavelet). 
The cleaned maps are combined in the Bayesian estimator and 
the output map of Compton parameters is derived. SExtractor is 
applied to the map of reconstructed Compton parameter to de- 
tect objects above a given threshold and compute their flux. The 
thresholds are based on the background or noise estimated by 
SExtractor. In order to compute the purity, different signal-to- 
noise ratios (ranging from 3 to 10) are used to compute the flux. 
The method assumes a power spectrum for the cluster population 
although this is not critical. 

The main advantage of the method resides in its robustness 
(almost no assumptions) and its ability to reconstruct both ex- 
tended and compact clusters. The main limitation is the rela- 
tively poor reconstruction of the total flux of the cluster as com- 
pared to matched filters. 

3.3. Internal Linear Combination Methods 

The internal linear combination (ILC) is a simple method for 
extracting one single component of interest out of multifre- 
quency observations. It has been widely used for CMB estima- 
tion on WMAP data (Bennett et al. 2003; Tegmarket al. 2003; 
Eriksen e t al. 2004; P ark et all 120071; IPelabrouille et al.1 120091: 
Kim et al.i.2009: Samal et a J2009I). A general descripti on of the 
method can be found in IPelabrouille & Cardosol (l2009h . 

The general idea behind the ILC is to form a linear combina- 
tion of all available observations which has unit response to the 
component of interest, while minimizing the total variance of the 
output map. This method assumes that all observations yi(p), for 



channel i and pixel p, can be written as the sum of one single 
template of interest scaled by some coefficient a,, and of unspec- 
ified contaminants which comprise noise and foregrounds, i.e. 

yiip) = ais(p) + ni(p) (5) 

where yiip) is the observed map in channel /, s(p) is the template 
of interest (here, the SZ map), and «,(/?) comprise the contribu- 
tion of both all astrophysical foregrounds (CMB, galactic emis- 
sion, point sources...) and of instrumental noise. This equation 
can be recast as: 

yip) = a sip) + nip) (6) 

To first order, linear combinations of the inputs of the form 
Tji^iyiip) guarantee unit response to the component of inter- 
est provided that the constraint 2, w,fl, =1 is satisfied (there 
are, however, restrictions and higher order effects, which are 
discussed in detail in the appendix of Pelabrouille et al. (2009!) 
and in |P ick et al. (201 0|)). It can be shown straightforwardly 
(lEriksenet al..,2004; .Pelabrouille & Cardosol l2009) that the lin- 
ear weights which minimize the variance of the output map are: 



a'R 'a 

where R is the empirical covariance matrix of the observations. 
What distinguishes the diff'erent ILC implementations is essen- 
tially the domains over which the above solution is implemented. 

3.3.1. ILC1 

The needlet ILC method works in two steps. First, an SZ map 
is produced by internal linear combination, with a n eedlet space 
implementation similar to that of Pelabrou ille et al.l (2009). The 
use of spherical needlets permits the ILC filter to adapt to lo- 
cal conditions in both direct (pixel) space and harmonic space. 
Input maps include all simulated Planck maps, as well as an 
external template of emission at 100 microns dSchlegel et al.l 
1998), which helps subtracting dust emission. The cluster cat- 
alogue is then obtained by matched filtering on small patches 
extracted from the needlet ILC SZ full-sky map, as described in 
iMelin etalJ(l2006h . 

3.3.2. ILC2 

The ILC2 approach relies on the same processing for photom- 
etry, cluster size, and signal-to-noise estimates as in ILCl, but 
in this case the detection of cluster candidates is made using 
SExtractor on a Wiener-filtered version of the ILCl map. 

While the extraction of the ILC map works on full sky maps, 
using the needlet framework to perform localized filtering, here 
again the detection of cluster candidates, and the estimation 
of size and flux, are performed on small patches (obtained by 
gnomonic projection). 

3.3.3. ILC3 

In this method, a filter in the harmonic domain is applied to con- 
struct a series of maps that are sensitive to the range of clus- 
ter scales. We used a Mexican-hat filter constructed from two 
Gaussians, one with 1/4 the width of the other A list of clus- 
ter candidates is compiled using a peak finding algorithm, which 
searches for enhanced signal levels in the individual map by fit- 
ting cluster model parameter. We employed the y6-model profile 
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with p - 2/3 convolved with the Planck beams. The catalogue 
produced then includes as parameters the cluster location, the 
flux, and the size estimate. The errors on these parameters are 
incorporated as given by the likelihoods of the fit. 

Additional improvements to the method can be achieved by 
using more optimal foreground estimators (but probably with 
slower convergence). This is useful especially for the fainter 
clusters, or those confused to a high degree by source contami- 
nation. 



3.3.4. ILC4 

The ILC4 method is a standard ILC in Fourier space performed 
on the square patches. It is implemented independently in an- 
nuli in wave number (modulus of the Fourier mode) by applying 
weights, according to Equation |7] this time in the Fourier do- 
main. The cluster detection is performed using SExtractor on the 
reconstructed SZ map. The flux estimation is performed on the 
original multifrequency maps (small patches) using the MMF at 
the position of the SExtractor detections. 

3.3.5. ILC5 

This algorithm is designed to work on local noisy multichannel 
maps in the wavelet domain. The representation of galaxy clus- 
ters in an appropriate biorthogonal wavelet basis is expected to 
be sparse compared to the contributions of other astrophysical 
components. This should ease the subsequent SZ separation us- 
ing the ILC method at each wavelet scale. The reverse wavelet 
transform is then applied to estimate local SZ-maps. The latter 
are convolved using a Gaussian beam of 5 arcmin FWHM to re- 
duce the noise prior to cluster detection using SExtractor with a 
threshold fixed classically to a multiple of the rms noise (signal- 
to-noise ratio ranging from 3 to 6). The brightest IR galaxies and 
radiosources are masked to reduce contamination. Finally, mul- 
tiple detections due to the overlap of local maps are removed. 
Multiscale ILC proved to be more efficient than regular ILC to 
remove large angular scale contamination on local map simu- 
lations. However, the cluster catalogue appears comparatively 
more contaminated which may be due to an imperfect cleaning 
of multiple detections. Also, the known point sources were only 
masked before detection with SExtractor Masking the observed 
sky maps earlier could improve the component separation using 
Multiscale ILC. Doing so may prevent a few very bright pix- 
els biasing the ILC parameter estimation but in turn raises the 
question of data interpolation across masked regions. The other 
implementations of the ILC do mask point sources before com- 
bining the maps and are thus not subject to this bias. 



coefficients in O : 

min|s,si^l|S<I»^ll + ^l|Y-AS||2 (8) 

where || . . . II2 is the L2 (Euclidean) norm. In lBobinetal] ( l2008l) . 
it was shown that sparsity enhances the diversity between the 
components thus improving the separation quality. The spec- 
tral signatures of CMB and SZ are assumed to be known. The 
spectral signature of the free-free component is approximately 
known up to a multiplicative constant (power law with fixed 
spectral index). 

Hence, GMCA furnishes a noisy SZ map in which we want 
to detect the SZ clusters. This is done in three steps; 

- wavelet denoising 

- SExtractor to extract the clusters from the noise-free SZ map, 
and finally 

- a maximum likehood to get the flux of the detected sources. 
4. Results 

We evaluated each extracted catalogue in terms of catalogue 
content and photometric recovery based on comparisons be- 
tween the extracted catalogues and the simulated input SZ clus- 
ter catalogue. For this purpose, we cut the input catalogue at 
y > 5 X 10 arcmin^, well below the theoretical Planck detec- 
tion limit (see below), and restrict ourselves to the high latitude 
sky at \b\ > 20 degs to reduce contamination by galactic fore- 
grounds. We then cross-match the candidate cluster in a given 
catalogue to a corresponding input cluster Each match results in 
a true detection, while candidates without a match are labeled 
as false detections. In a second step, we compare the extracted 
properties, namely SZ Compton parameter and size, of the true 
detections to the input cluster properties. 

Angular proximity was the only association criterion used 
for the matching. Specifically, we matched an extracted cluster 
to an input cluster if their separation on the sky 6 < Omax = fiOv), 
a function of the true angular virial radius of the (input) cluster 
The function /(0,.) varied over three domains: f(6y) - 5 arcmins 
for 0,, < 5 arcmins; f(0v) = dv for 5 arcmins < 0,, < 20 arcmins; 
and f{Oy) = 20 arcmins for 0^ > 20 arcmins. 

We first focus on the catalogue completeness and purity, 
both of which we define immediately below. We then test the 
accuracy of the recovered flux and size estimates, as well as 
each algorithm's ability to internally evaluate the uncertainties 
on these photometric quantities. Since many fewer codes ran on 
the Challenge v2, we focus mainly on Challenge vL We include 
results for those codes that did run on Challenge 2 to gauge the 
influence of the underlying cluster model used for the simula- 
tion. 



3.4. GMCA 

Generalized morphological component analysis is a blind source 
separation method devised for separating sources from instanta- 
neous linear mixtures using the model given by; Y - AS + n. 
The components S are assumed to be sparsely represented (i.e. 
have a few significant samples in a specific basis) in a so-called 
sparse representation O (typically wavelets). Assuming that the 
components have a sparse representation in the wavelet domain 
is equivalent to assuming that most components have a certain 
spatial regularity. These components and their spectral signa- 
tures are then recovered by minimizing the number of significant 



4.1. Catalogue Content 

A useful global diagnostic is the cur ve of y/e/t/ v ersus global 
purity for a given catalogue (see e.g. iPires et aT 1 (120061) ). The 
former is simply the total number of clusters detected and the 
latter we define as 1 - F^, where Fg is the global contamination 
rate; 

total number of false detections 
^ total number of detections 

The yield curve is parametrized by the effective detection thresh- 
old of the catalogue construction algorithm. It is a global diag- 
nostic because it gauges the total content of a catalogue, rather 
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Fig. 1. For SZ Challenge vl: Yield as a function of global purity. The right handside panel is a zoom on the high-purity region. 
Each curve is parameterized by the detection threshold of the corresponding algorithm. As discussed in the text, the overall value of 
the yields should be considered with caution, due to remaining modeling uncertainties (see text). We focus instead on relative yield 
between algorithms as a measure of performance. 




Fig. 2. For SZ Challenge v2: Yield as a function of global purity. The same comment applies concerning the overall yield values; in 
particular, the cluster model changed significantly between the versions v 1 and v2 of the Challenge resulting in lower overall yields 
here. Fewer codes participated in the Challenge v2 (see text). 



than its content as a function of flux or other measurable quanti- 
ties. 

Figure [T] compares the yield curves of outputs of all the 
algorithms in the SZ Challenge vl, and Figure |2] those for 
the SZ Challenge v2. Increasing the detection threshold moves 
a catalogue along its curve to higher purity and lower yield. 
Algorithms increase in performance towards the upper right- 
hand corner, i.e., both high yield and high purity. 

As to be expected, algorithms that locally estimate the noise 
(both instrumental and astrophysical), i.e. on local patches of 
a few square degrees, perform much better than those that em- 
ploy a global noise estimate, such as MMF4 and ILC3. For those 
methods with local noise estimation, we note that their effec- 
tive survey depth appears to anticorrelate with the instrumental 
noise, indicating that astrophysical confusion is effectively re- 



moved. This can be seen in Figure |3] which compares the den- 
sity of detected SZ sources (top panel) to the pixel hit count (bot- 
tom panel). The result illustrated with one single method, ILC2 
run on the SZ challenge vl, holds for the other algorithms. The 
cluster detection limit appears to be primarily modulated by the 
instrumental noise at high Galactic latitude, as opposed to fore- 
ground emission. 

Less expected, perhaps, is the fact that all algorithms tend to 
miss nearby clusters. These are extended objects, and although 
they have large total SZ flux, these clusters are "resolved-out" 
- missed because of their low surface brightness. This is an ex- 
treme example of resolution effects expected in the case of SZ 
detection in relatively low resolution experiments like Planck. It 
is not related to the foreground removing efficiency since the ef- 
fect can also be mimicked in simulations including only instru- 
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Fig. 3. Illustrated for ILC2 on the SZ Challenge vl, the detection density (top panel) is compared to the pixel hit count for the map 
at 143 GHz (bottom panel). The noise in the simulated Planck maps scales as 1 / yjhits. At high Galactic latitudes, the detection 
density clearly anticorrelates with map noise. Both maps are smoothed on a scale of 20 degrees, and the Galactic plane (\b\ < 20 
deg) is masked. 



mental white noise. Apart from the few clusters that are fully 
resolved, many have angular sizes comparable to the ef fective 
beam , and this leads to a non -trivial selection function (IWhitd 
I200I iMelin etani2005[ l2006h . 

We emphasize that the numerical values of the yields de- 
pend on the cosmological model, on the foreground model and 
on the cluster model used in the simulations. They must be con- 



sidered with caution because of the inherent modeling uncer- 
tainties. As for the foreground model, the templates used to 
model Galactic components in the PSM were chosen so that 
they are reasonably representative of the complexity of the dif- 
fuse galactic emission. Thanks to man y new observa t ions in 
particular in the IR and submm doma in dLagache et al. (2007), 
Viero et al.. (,2009i) . iHall et all (l20Tol) . lAmblard et al., (.20 111) . 
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iPlanck Collaboration et alJ (|201 id) '), the models of point sources 
have evolved very much between the beginning of the SZ chal- 
lenge and the publication of these results. These updates were 
not taken into account in the PSM when the study was per- 
formed. Moreover, the cluster model in challenge vl was based 
on the isothermal yS-model, while v2 employed a modified NFW 
pressure profile favoured by X-ray determinations of the gas 
pressure (Arnaud et al. 2010) with a normalization of the Y - M 
relation lower by ~ 15% than in vl. Finally, erg changed from 
0.85 in challenge vl to 0.796 in v2 which strongly influences 
the total cluster yield. 

The more peaked profile actually improves detection ef- 
ficiency, while the lower normalization reduces the predicted 
yield. Along with the much lower value of era, the net result is 
that the yields in challenge v2 are noticeably lower than in vl 
as seen in Figs.[T]and|2] We thus only discuss, in this study, the 
relative yields of the codes as a gauge of performance treating 
the absolute value of the yield with caution. 

Focusing on the relative merit of the algorithms, we see that 
Figures [1] and |2] display large dispersion in the yield at a given 
purity. This reflects of course the intrinsic performance of the 
algorithms, but also for the detection methods that share simi- 
lar underlying algorithms, e.g. MMF and ILC, the dispersion in 
the yield reflects the differences in implementations (e.g. noise 
estimation, de-blending, etc). 

Somewhat deceptively, these yield variations correspond to 
only minor differences in detection threshold, as illustrated for 
the SZ Challenge vl in Figure |4] This figure traces the curves 
above which 90% and 10% (lower set and upper set respectively) 
of the clusters detected by each method lie in the true Y - true 
0y plane. As already mentioned, many clusters are marginally 
resolved (sizes at least comparable to the beam), which means 
that detection efficiency depends not only on flux, but also on 
siz^ 

The algorithms all have similar curves in this plane. This 
means that the differences in yield are due to only small vari- 
ations of the selection curve since completeness is expected to 
be monotonic. The black points represent a random sample of 
1/4 of the input clusters and show where the bulk of the cata- 
logue lies. These small variations have important consequences 
for cosmological interpretation of the counts and hence must be 
properly quantified. 

We compute for all the methods and in the SZ challenges 
vl and v2 the completeness defined as the ratio [true detections 
( recovered clusters)/simulated clusters] over bins of true (simu- 
lated) flux Y. We find it varies from 80 to 98% at Yn^ = 10"^ 
arcmin^ for the direct methods (based on frequency maps). The 
completness is of order of 80% at the same Fiim for the indirect 
methods (based on detections in SZ maps). We note a slight in- 
crease of the completeness from the SZ Challenge vl to v2. We 
also estimate the contamination of the output catalogues defined 
as the ratio [false detections/total detections]. This is evaluated 
as as a function of recovered flux. The average contamination 
of the output catalogues ranges between 6 and 13% both in the 
case of the challenge 1 and 2. However, the purity with respect 
to the Y bins differ significantly from method to method. As a 
general trend, the lowest Y bin, i.e. the smallet recovered fluxes 
centered around lO^^* arcmin^, is more contaminated (~75% on 
average) in the case of the indirect method than in the case of 
direct methods (~50% on average). 



* This was shown on real data in lPlanck Collaboration et aP (1201 Ibl) 



4.2. Photometric and Astrometric Recovery 

Cluster characterization is a separate issue from detection. It in- 
volves determination of angular positions as well as photometry. 
Since Planck will marginally resolve many clusters, photometry 
here means both flux Free and characteristic size measurements. 
Moreover, each method should provide an estimate of the errors 
on these quantities for each object in the catalog. 

We illustrate in Fig. |5] a scatter diagram of positional off'- 
set for MMF3 as a function of true cluster size, 0y. On average, 
all the algorithms perform similarly and recover cluster position 
to ~ 2 arcmins with a large scatter In addition, we see that it 
is more difficult to accurately determine the position of intrin- 
sically extended clusters, as shown by the fact that the cloud of 
points is elongated and inclined. 

Concerning photometry, we show in Figure|6]and|7]the mean 
recovered SZ flux. Free, normalized to the true (simulated) flux, 
F, as a function of the latter Only true detections are used in this 
comparison. We illustrate our results for a subset of methods 
namely MMFl, MMF2, MMF3, ILCl, 1LC2, 1LC3. Methods 
that filter the maps such as ILC4 and GMCA, or that use the 
SExtractor photometry such as BNP and 1LC5 exhibit significant 
bias in flux recovery at the bright end. At the faint end, we see the 
appearance of Malmquist bias as an upturn in the measured flux. 
The importance of this bias varies from algorithm to algorithm. 

Figure |7] gives the dispersion in the recovered flux cry^^^ as a 
function of true F (once again, only involving true detections). 
Here, we see that some codes perform significantly better than 
others. Those that adjust an SZ profile to each cluster outper- 
form by a large margin those that do not. Photometry based 
on SExtractor, for example, fares much worse than the MMF 
codes. Even among the best performing codes, however, the in- 
trinsic photometric dispersion is of order 30%. This is impor- 
tant, because we expect SZ flux to tightly correlate with clus- 
ter mass, with a scatte r as low as ~ 10% as indicated by both 
numerical simu lations I da Silva et"al] l2004b; Motl et al. 2005S 
Kravtsov et al. 2006) and recent X-ray data (Arnaud et al. 2007t 
Nagai et al. 2007). The SZ flux hence should offer a good mass 
proxy. What we see from this figure, however, is that the obser- 
vational scatter will dominate the intrinsic scatter of this mass 
proxy and needs to be properly accounted for in the cosmologi- 
cal analyses. 

We have attributed the origin of the photometric scatter to 
difficulty in determining cluster size. Although methods adjust- 
ing a profile to the SZ are able to estimate the size of many clus- 
ters, they do so with significant dispersion. Furthermore, this is- 
sue arises specifically for Planck-like resolution because a large 
number of clusters are only marginally resolved. Imposing the 
cluster size, for example from external data, such as X-ray or op- 
tical observations or higher resolution SZ measurements, would 
significantly reduce the observational scatter. 

5. Discussion and conclusion 

In the present study, we compare different codes and algorithms 
to detect SZ galaxy clusters from multi-wavelength experiments 
using Plank's instrumental characteristics. These methods may 
be usefully divided into direct methods (four matched-filter ap- 
proaches and PowellSnakes) using individual frequency maps, 
and indirect methods (five ILC methods, GMCA and BNP) that 
first construct an SZ map in which they subsequently search for 
clusters. 

As already emphasized, the global yield values of all meth- 
ods must be considered with caution because of inherent model- 
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Fig. 4. Selection curves for each algorithm in the true Y - true 6^ plane for SZ Challenge vl. The lower set of curves indicate the 
90th percentiles, i.e., the curve above which lie 90% of the detected clusters; the upper set corresponds to the 10th percentile. The 
color codes are as in Figure [1] We see that the Planck selection function depends not just on flux, but also on cluster angular size. 
Many clusters are at least marginally resolved by Planck, leading to these size effects in the selection fuction. The dashed lines show 
contours of fixed mass and redshift, as indicated, while the cloud of points shows the distribution of the input catalogue(in fact a 
subsample of 1 in 4 randomly selected input clusters). We see that small variations in selection curves generate significant yield 
changes. 



ing uncertainties of the sky simulations and cluster models used, 
and to the underlying cosmological model. Therefore, we focus 
on relative yields as a gauge of performance of the algorithms. It 
is worth noting that results of a direct or indirect method signif- 
icantly vary (within factors of as much as three) with the details 
of their implementation, with clear impact on the survey yield 
as demonstrated in Figures[T]and|2l Using the PSM simulations 
and including th noise as described in 12.11 we would expect of 
order of 1000-2000 clusters at a purity of ~ 90% with \b\ > 20 
deg. This number depends on the extraction method used and 
may vary with a more detailed modeling of the sky. The cluster 
yield can be increased by accepting a higher contamination rate 
and calling for extensive follow-up to eliminate false detections 
a posteriori. 

The indirect methods seem to offer greater opportunity for 
optimization with a larger number of tuning parameters. They 
are also less model dependent for the SZ map construction 
and the cluster detection. Although, they can be coupled with 
matched filters for the SZ flux measurement. In turn, the direct 
methods are linear, easy to implement and robust. One of their 
advantages relies in the fact that they can be optimized to detect 
objects of a given shape (SZ profile) and and a given spectral 
energy distribution - SED - (SZ spectrum). This characteristic 
of the direct method is particularly important for Planck-like 
multi-frequency surveys with moderate resolution. Indeed, due 
to lack of resolution, spurious sources (galactic features or point 
sources) may be detected by the spatial filter and in that case the 
frequency coverage and the spectral matching is the best strat- 
egy to monitor the false detections. Due to their robustness and 
easy implementation direct methods are more adapted to run in 



pipelines^ The situation is quite different for high resolution, 
arcminute-scale, SZ experiments such as ACT and SPT where 
the filtering of one unique low frequency map (where the SZ 
signal is negative) is sufficiently efficient to unable cluster detec- 
tion. In these cases though, extended clusters are not well recov- 
ered as they suffer more from the CMB contamination and thus 
from the filtering. 

The comparison of different codes and cluster detection 
methods exhibits selection effects and catalogue uncertainties, 
neither of which depend, for example, on the actual cluster 
physics model. This is shown in the selection curves Fig. |4l 
In a Planck-like case, with moderate resolution, clusters do not 
appear as point sources, but many are resolved or have sizes 
comparable to the effective SZ beam. In view of this, the use 
of an adapted spatial filter to optimally model the SZ profile 
provides a significant improvement in the detection yield and in 
the photometry. Clusters being partially resolved leads to non- 
trivial selection criteria that depend both on flux and true an- 
gular size, as demonstrated by the fact that the curves in Fig. |4] 
are not horizontal lines. In that respect, the use of X-ray informa- 
l :ion from the ROSAT All-Sky Survey (RASS) clus ter catalogues 
jBohringer etal]|2000ll2004t lPiffaretti et al."201 1) or from opti- 
cal catalogues, e.g. in the SDSS area (Koester et al. 2007), will 
be of particular value as they will give us a handle on both flux 
and size of the clusters detected by Planck and, even more im- 
portantly, understand better the completeness in studying those 



' MMFl, MM2, MMF3 and PwS were all implemented in Planck's 
pipelines. Furthermore, MMF3 was used to extract the Planck clusters 
published in the ESZ sample. 
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Fig. 5. Positional accuracy of the recovered clusters illustrated for MMF3 in the SZ challenge v2. 




Fig. 6. Flux recovery bias. The Figure shows for a subset of methods the average recovered Free, normalized to the true input Y, as a 
function of Y. At the bright end, most codes extract an unbiased estimate of cluster flux, while the expected Malmquist bias appears 
at the faint-end, just below Y ~ 2x 10"^ arcmin^. 



which are missed dChamballu et alfcOlOl) . As already stated, the 
exact position of each selection curve in Fig. |4] depends on the 
algorithm and small variations in the position of this curve pro- 
duce significant changes in catalogue yield. Most of the differ- 
ences bewteen catalogues occur at small flux and size, where 
the bulk of the cluster population resides. These objects are for 
the majority low mass, intermediate redshift clusters. Moreover, 
detection becomes progressively less efficient for large objects; 
this is intrinsic and hence shared by all algorithms. 



Concerning individual cluster measurements we find that the 
astrometry is recovered to ~ 2 arcmins on average and photome- 
try to ~ 30% for the best-behaving algorithms. The positional er- 
ror is not a problem for X-ray/optical follow-ups because Planck 
is expected to detect massive clusters which will be easy to find 
in the XMM/Chandra/1 to 4meter class optical telescope fields 
of view. The photometric error in our Planck-like case is domi- 
nated by the difficulty in accurately determining the cluster ex- 
tent/size. This is another consequence of the fact that many clus- 
ters are marginally resolved by Planck: large enough that their 
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Fig. 7. Flux recovery uncertainty for the subset of methods given in Figure |6] The Figure shows the dispersion in measured flux 
about the true input Y flux as a function of the true input flux. The best algorithms are affected by ~ 30% dispersion. 



angular extent matters, but small enough that we have difficulty 
fixing their true size. One way of reducing the photometric er- 
ror is thus using external constraints on cluster size. One again 
ancillary data from RASS or optical cluster catalogues will help 
in this regard, at least at low redshift (z < 0.3 - 0.5) ; at higher 
redshift, we will rely on follow-up observations if we want to 
reduce photometric uncertainties. 

The comparison of an ensemble of cluster extraction meth- 
ods in the case of a multi-frequency moderate resolution experi- 
ment shows that the optimization of the cluster detection in terms 
of yield and purity, but also in terms of positional accuracy and 
photometry, is very sensitive to the implementation of the code. 
The global or local treatment of the noise estimate or the clean- 
ing from point sources are the two main causes of difference. 
However and most importantly, the use of as realistic as possi- 
ble SZ profile (as opposed to model independent profile) to filter 
out the signal or to measure the fluxes is a key aspect of cluster 
detection techniques in our context. In that respect, using exter- 
nal information from SZ observations or from other wavelength 
will significantly help in improving the measurement of the clus- 
ter properties and in turn optimize the catalogue yields and their 
selection function. 
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Appendix A: Comparison plots for Cliallenge SZ v1 

The plots are obtained at purity=0.9 (except when the catalog 
does not permit this value to be reached). 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



10" 

o , 

10^ 

10° 



Positional accuracy 



10" 









Simulated (13747 clusters) : 

Recovered (2272 clusters) ' 







10" 



10" 



10" 



10"' 10° 



Y simulated [arcmin^] 

Completeness 




lO""* 10"^ 10"^ 10"' 
Y simulated [arcmin^] 

SZ flux 




"D 
0) 

o 

E 
'in 

>- 



(U 
L_ 

> 

o 
u 

0) 



o 



10"" 10"^ 
Y simulated [orcmin^] 

Mean &c variance 




10"^ 10""* 10"^ 10"^ 10"' 10° 
Y simulated [arcmin^] 




10 100 1000 

e„ simulated [arcmin] 

Contannination 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 



1.0 
0.8 

0.6 

0.4 

0.2 
0.0 



Totol=9.7% 




10" 



^ 1000 
c 

E 
o 

^ 100 

T3 
(U 

i_ 
0} 

o 



10"" 10" 



10" 



10" 



10° 



Y recovered [arcmin^] 

SZ viriol radius 



1 




1 10 100 1000 

0„ simulated [arcmin] 



^— ' 
o 

D 

.E 

'm 

OS* 

\ 
x> 

0) 
0) 

> 
o 
o 



cn 
o 



Mean &c variance 



1.0 

0.5 

0.0 

-0.5 
-1.0 



1 10 100 1000 

9^ simulated [arcmin] 



Table A.l. MMFl 



16 



J.-B. Melin et al.: A Comparison of Algorithms for the Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



10" 
10^ 

o , 

10^ 

€ 10^ 

10^ 
10° 



Positional accuracy 



10" 



Simulated (13747 clusters) 
Recovered (1891 clusters) 




10" 



10" 



10" 



Y simulated [arcmin^] 

Completeness 



lO""* 10"^ 10"^ 10"' 
Y simulated [arcmin^] 

SZ flux 



10"'' 10"^ 
Y simulated [orcmin^] 



10"' 10° 






10 100 
e„ simulated [arcmin] 

Contonnination 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 



Y recovered [arcmin^] 

SZ viriol radius 



^ 1000 



c 

E 
o 



^ 100 



T3 
(U 

0} 

> 
o 
o 



1000 





1 10 100 1000 

0„ simulated [arcmin] 



"D 
(U 

o 

E 
'in 

>- 



(U 
L_ 

> 

o 
u 

0) 



o 



Mean &c variance 




10"^ lO"'' 10"^ 10"^ 10"' 10° 
Y simulated [arcmin^] 



-I—' 
o 

D 

.E 

'm 

ccT 

\ 
x> 

0) 
0) 

> 
o 
o 



cn 
o 



Mean &c variance 




1 10 100 1000 

9^ simulated [arcmin] 



Table A.2. MMF2 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



10" 

o , 

10^ 

10° 



10" 



10" 
10"' 



^ 10-^ 

"O 



10" 



10" 



>- 10" 



10" 



10" 



10" 



Y simulated [arcmin^] 

Completeness 



10" 



lO""* 10"^ 10"^ 10"' 
Y simulated [arcmin^] 

SZ flux 









Simulated (13747 clusters) : 

Recovered (1701 clusters) 







10"' 10° 






'y-. 

y ■■ 




y 
<^ 










' ' 





10"^ 10"' 



10"" 10"^ 
Y simulated [orcmin^] 



10" 



Positional occurocy 




10 100 
e„ simulated [arcmin] 

ContGnninotion 



in 
c 
o 

o 

(D 

0) 
Q 
\ 

(D 

in 



Y recovered [arcmin^] 

SZ viriol radius 



^ 1000 



c 

E 
o 



^ 100 



T3 
(U 

> 

o 
o 



1000 





1 10 100 1000 

0„ simulated [arcmin] 



"D 
(U 

o 

E 

>- 



"D 
(U 
L_ 

> 

o 
u 

0) 



o 



Mean &: variance 




10 



10""* lO"-" 10"^ 10" 

Y simulated [arcmin^] 



10° 



-I—' 

D 

.E 

OS* 

\ 

x> 

0) 
0) 

> 
o 
o 



cn 
o 



Mean &: variance 



1.0 

0.5 

0.0 

-0.5 
-1.0 





T ' ' ' — 1 ' ' ' — 




,||j4il^i{ 



1 10 100 1000 

0„ simulated [arcmin] 



Table A.3. MMF3 



18 



J.-B. Melin et al.: A Comparison of Algorithms for the Construction of SZ Cluster Catalogues 



Simulated &c recovered clusters 



10'' 
10^ 

;o 10^ 

? 10^ 

10' 
10° 



10" 



Simuloted {13747 clusters) 
Recovered (863 clusters) 




lO""* 10" 



10" 



10" 



Y simulated [arcmin ] 



10° 



c 
E 

u 



u 
o 

D 

o 
u 

D 

"5 

C 

_o 

"Si 

o 



Positional accuracy 



10.00 



1.00 r 



0.10 



0.01 




10 100 1000 

0, simulated [arcmin] 



Completeness 




Y 



10 " 10 ' 10 
simulated [orcmin^] 



10" 



Table A.4. MMF4 (photometry has not been computed for this 
algorithm. 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



10" 

o , 

10^ 

10° 



Positional accuracy 



10" 



E 



"O 

(D 
> 

o 
o 



>- 



10" 



10" 



10" 



Y simulated [arcmin^] 

Completeness 



lO""* 10"^ 10"^ 10"' 
Y simulated [arcmin^] 

SZ flux 









Simulated (13747 clusters) : 

Recovered (1477 clusters) 







10"' 10° 





10"" 10"^ 
Y simulated [orcmin^] 




10 100 1000 

e„ simulated [arcmin] 

Contonnination 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 




Y recovered [arcmin^] 

SZ viriol radius 



^ 1000 



c 

E 
o 



^ 100 



T3 
(U 

0} 

> 
o 
o 




1 10 100 1000 

0„ simulated [arcmin] 



"D 
(U 

o 

E 
'in 

>- 



(U 
L_ 

> 

o 
u 

0) 



o 



Mean &c variance 




10"^ 10""* 10"^ 10"^ 10"' 10° 
Y simulated [arcmin^] 



-I—' 
o 

D 

E 

'm 

ccT 

\ 
x> 

> 

o 
o 



cn 
o 



Mean &c variance 



1.0 

0.5 

0.0 

-0.5 
-1.0 



1 1 1 1 1 r 


— 1 — 1 1 — 1 — 1 — 

: 


i , . , . 1 . . . . 





1 10 100 1000 

9^ simulated [arcmin] 



Table A.5. PS 



20 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



10" 

o , 

10^ 

10° 



Positional accuracy 



10" 



10" 



10" 



10" 



Y simulated [arcmin^] 

Completeness 



lO""* 10"^ 10"^ 10"' 
Y simulated [arcmin^] 

SZ flux 









Simulated (13747 clusters) : 

Recovered (1853 clusters) 







10"' 10° 





10"" 10"^ 
Y simulated [orcmin^] 




10 100 1000 

e„ simulated [arcmin] 

Contonnination 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 




Y recovered [arcmin^] 

SZ viriol radius 



^ 1000 



c 

E 
o 



^ 100 



> 

o 
o 




1 10 100 1000 

0„ simulated [arcmin] 



"D 
(U 

o 

E 
'in 

>- 



(U 
L_ 

> 

o 
u 

0) 



o 



Mean &c variance 




10"^ 10""* 10"^ 10"^ 10"' 10° 
Y simulated [arcmin^] 



-I—' 
o 

D 

.E 

'm 

OS* 

\ 
x> 

0) 
0) 

> 
o 
o 



cn 
o 



Mean &c variance 




1 10 100 1000 

9^ simulated [arcmin] 



Table A.6. ILCl 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



>- 

O 

T3 



10" 
10^ 
10'' 
10^ 
10^ 

10' 
10° 



Positional accuracy 



10" 









Simulated (13747 clusters) : 

Recovered (2143 clusters) 







10" 



10" 



10" 



10"' 10° 



Y simulated [arcmin^] 

Completeness 




Y simulated [arcmin^] 

SZ flux 




"D 
(U 

o 

E 
'in 

>- 



(U 
L_ 

> 

o 
u 

0) 



o 



10"^ 10"" 10"^ 10"^ 10" 



Y simulated [orcmin^] 

Mean &c variance 



10" 




10"^ 10"" 10"^ 10"^ 10"' 10° 

Y simulated [arcmin^] 




10 100 
e„ simulated [arcmin] 

Contannination 



IOC 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 




Y recovered [arcmin^] 

SZ viriol radius 



^ 1000 
c 

E 
o 

^ 100 



T3 
(U 

0} 

> 
o 
o 



10 



1 




-I—' 
o 

D 

.E 

'm 

OS* 

\ 
x> 

0) 
0) 

> 
o 
o 



cn 
o 



1 10 100 1000 

0„ simulated [arcmin] 

Mean &c variance 




10 100 1000 

9^ simulated [arcmin] 



Table A.7. ILC2 



22 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted & recovered clusters " 



10" 

o , 

10^ 

10° 



10" 



10" 



10" 



10" 



Y simulated [arcmin^] 

Completeness 



lO"'' 10"^ 10"^ 10"' 
Y simulated [arcmin^] 









Simulated (13747 clusters) : 

Recovered (942 clusters) 




J ' U 1 



10"' 10° 




Positional occurocy 




10 100 1000 

e„ simulated [arcmin] 

ContGnninotion 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 




10"^ 10"" 10"^ 10"^ 10"' 10° 
Y recovered [arcmin^] 



SZ flux 



SZ viriol radius 




10"" 10"^ 

Y simulated [arcmin^] 



^ 1000 
c 

E 
o 

^ 100 



T3 
(U 

0} 
> 

o 
o 



10 



1 





M 




y 

y 

y 




y 






^^^^ 




y 

y ' 





1 10 100 1000 

0„ simulated [arcmin] 



"D 
(U 

o 
E 



"D 
(U 
L_ 

> 

o 
u 

0) 



o 



Mean &: variance 




Y simulated [arcmin ] 



-I—' 

D 

.E 

OS* 

\ 

x> 

0) 
0) 

> 
o 
o 



cn 
o 



Mean &: variance 




1 10 100 1000 

0„ simulated [arcmin] 



Table A.8. ILC3 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



>- 

O 

T3 



10" 
10^ 
10'' 
10^ 
10^ 

10' 
10° 



10" 



Simulated (13747 clusters) 
Recovered (1789 clusters) 




10" 



10" 



10" 



10"' 10° 



Y simulated [arcmin^] 

Completeness 




Y simulated [arcmin^] 

SZ flux 




"D 
0) 

o 

E 
'in 

>- 



(U 
L_ 

> 

o 
u 

0) 



o 



10"^ 10"" 10"^ 10"^ 10" 



Y simulated [orcmin^] 

Mean &c variance 



10" 




10"^ 10"" 10"^ 10"^ 10"' 10° 

Y simulated [arcmin^] 



E 
o 



o 

D 

□ 
O 

o 
o 

o 
c 
_o 

'm 
o 



Positional occurocy 



10.00 

1.00 

0.10^ 
0.01 




1 10 100 1000 

e„ simulated [arcmin] 

Contannination 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 




10"^ 10"" 10"^ 10"^ 10"' 10° 
Y recovered [arcmin^] 



SZ viriol radius 



^ 1000 
c 

E 
o 

^ 100 



T3 
(U 

0} 

> 
o 
o 



10 



1 








r <> * 










o 






' 





(D 
-I—' 
O 
D 

E 

OS* 

\ 

x> 

0) 
0) 

> 
o 
o 



cn 
o 



1 10 100 1000 

0„ simulated [arcmin] 

Mean &c variance 




1 10 100 1000 

0„ simulated [arcmin] 



Table A.9. ILC4 



24 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



10" 

o , 

10^ 

10° 



10" 



10" 



10" 



10" 



Y simulated [arcmin^] 

Completeness 



lO""* 10"^ 10"^ 10"' 
Y simulated [arcmin^] 









Simulated (13747 clusters) : 

Recovered (1129 clusters) 




[-'"^ — h -1 1 



10"' 10° 




Positional occurocy 




10 100 1000 

e„ simulated [arcmin] 

Contonnination 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 




Y recovered [arcmin^] 



10" 
10"' 



^ 10"^ 

"O 



10" 



10" 



>- 10" 



10" 



SZ flux 





'y-. 

y ■■ 




y 















y 




y 

y 

y 




y 

' ' 





10"^ 10"' 



10"" 10"^ 
Y simulated [arcmin^] 



10" 



"D 
(U 

O 

E 
'in 

>- 



(U 
L_ 

> 

o 
u 

0) 



o 



Mean &c variance 



1.0 

0.5 

0.0 

-0.5 
-1.0 




10"^ 10""* 10"^ 10"^ 10"' 10° 
Y simulated [arcmin^] 



Table A.IO. ILC5 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



10" 

o , 

10^ 

10° 



10" 



10" 



10" 



10" 



Y simulated [arcmin^] 

Completeness 



lO""* 10"^ 10"^ 10"' 
Y simulated [arcmin^] 

SZ flux 









Simulated (13747 clusters) : 

Recovered (1278 clusters) 




H"^ nil—, ^ 



10"' 10° 





10"" 10"^ 
Y simulated [orcmin^] 



Positional occurocy 




10 100 1000 

e„ simulated [arcmin] 

ContGnninotion 



m 
c 
o 

o 

(D 
-t—f 

0) 
Q 
\ 

(D 

m 




Y recovered [arcmin^] 

SZ viriol radius 



^ 1000 
c 

E 
o 

^ 100 



T3 
(U 

0} 

> 
o 
o 



10 



1 




W> <nMM»0O 00 



1 10 100 1000 

0„ simulated [arcmin] 



"D 
(U 

o 

E 
'in 

>- 



(U 
L_ 

> 

o 
u 

0) 



o 



Mean &c variance 




10"^ 10""* 10"^ 10"^ 10"' 10° 
Y simulated [arcmin^] 



-I—' 
o 

D 

E 

'm 

ccT 

\ 
x> 

0) 
0) 

> 
o 
o 



cn 
o 



Mean &c variance 




1 10 100 1000 

9^ simulated [arcmin] 



Table A.ll. GMCA 



26 



J.-B. Melin et al.: A Comparison of Algorithms for tiie Construction of SZ Cluster Catalogues 



Simuloted &c recovered clusters " 



10" 

o , 

10^ 

10° 



10" 



10" 
10"' 



^ 10-^ 

"O 



10" 



10" 



>- 10" 



10" 



10" 



10" 



Y simulated [arcmin^] 

Completeness 



10" 



lO""* 10"^ 10"^ 10"' 
Y simulated [arcmin^] 

SZ flux 



10"^ 10"' 



10"" 10"^ 
Y simulated [orcmin^] 









Simulated (13747 clusters) : 

Recovered (1678 clusters) 







10"' 10° 




















y 

y 

y 










y 










y 










y 










fig^|H^i^o coo 






























n 










o 






y 

/ 










y 










/ 




















y 

y 



















10" 



Positional occurocy 




10 100 1000 

e„ simulated [arcmin] 

ContGnninotion 



in 
c 
o 

o 

(D 

0) 
Q 
\ 

(D 



1.0 

0.8 

0.6 

0.4 

0.2 
0.0 



Tot0l=13.4% 



1 



10"^ 10"" 10"^ 10"^ 10"' 
Y recovered [arcmin^] 

SZ viriol radius 



^ 1000 
c 

E 
o 

^ 100 



0) 
> 

o 
o 



10 



1 



10° 




1 10 100 1000 

0„ simulated [arcmin] 



"D 
(U 

o 
E 



"D 
(U 
L_ 

> 

o 
u 

0) 



o 



Mean &: variance 




10"^ 10"" 10"^ 10"^ 10"' 10° 
Y simulated [arcmin^] 



D 

.E 

OS* 

\ 

x> 

> 

o 
o 



cn 
o 



Mean &: variance 



1.0 

0.5 

0.0 

-0.5 
-1.0 



1 10 100 1000 

0„ simulated [arcmin] 



Table A.12. BNP 



J.-B. Melin et al.: A Comparison of Algorithms for the Construction of SZ Cluster Catalogues 27 

Appendix B: Comparison plots for Challenge SZ v2 
(Run 2) 

The plots are obtained at purity=0.9 (except when the catalog 
does not permit this value to be reached). 
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